model = 3;  % select pricing model (2-Calvo; 3-Golosov-Lucas; 4-Woodford, 5-linear hazard, 6-asymmetric linear hazard, 7 - uniform random menu cost plus)
mymoment = 'hazard';  %'stdpcdist 'hazard_USvol' 'hazard'; or 'kurtosis'  %'hazard'; %The moment to match
switch_distribution = 'gaussian';   %'gaussian' or 'laplace' or 'mixed'
nump = 501;   % Number of price grid points


estim = 0;  % estimate model parameters?
estim_dyn = 0; %estimate dynamic parameters?
doss  = 1;  % do steady state?
dodyn = 0;  % do dynamics?

estim_type = 'quick'; %'quick' requires fairly good starting value or 'robust'
switch_sim = 'irf';   %impulse response calculation 'irf' or simulation 'sim' 
switch_shock = 'm';   %shock type: 'A' or 'm'

freq_target = 0.1094;
size_target = 0.0896;
kurt_target = 2.71;
N_target    = 1/3; 
std_pi_yoy_target = 0.0095;   
std_pi_target = 0.02295;   
autocorr_pi_yoy_target = 0.9684;  %unused 

if ~exist('deviation_data')  %do it only the first time
            read_data;
            %recalculate frequency and size targets; exclude first and last
            %bins
            range_data = 2:length(density_data)-1;
            freq_target = density_data(range_data)'*freq_data(range_data)/sum(density_data(range_data));
            size_target = (density_data(range_data).*freq_data(range_data))'*abs(deviation_data(range_data))/sum((density_data(range_data).*freq_data(range_data)));
            mean_density_data = density_data(range_data)'*deviation_data(range_data)/sum(density_data(range_data));
            std_density_data = sqrt(density_data(range_data)'*(deviation_data(range_data)-mean_density_data).^2/sum(density_data(range_data)));
            kurt_density_data = density_data(range_data)'*(deviation_data(range_data)-mean_density_data).^4/(sum(density_data(range_data)).*std_density_data.^4);
end

noise = 0.0003;  %0.0003

alpha0 = NaN; menucost0 = NaN; omega0 = NaN; alphapos0=NaN;
lam0 = 1;
std_m0 = 0.003321603519177;
phiR0 = 0.61;  % Money growth shock
phiA0 = 0.95;   % TFP shock
std_A0  = 0.0013;

lbar=freq_target; % initialize parameters

% set initial GUESSES
switch model
    case 2  % Calvo
            std_eA0 = 0.041196011;  
            chi0 = 0.68513627;             
            wbar0 = 0.65656638;  
            chi0_min = 0.8;
            chi0_max = 1.1;
            wbar0_min = 0.3;
            wbar0_max = 0.4;
    case 3  % Golosov-Lucas
            menucost0 = 0.018751384;  %0.012479025;  %theta3:   
            std_eA0 = 0.029985967;    %0.029944736; %theta3:   
            chi0 = 0.66790114;   %0.50061324;  %%theta3:   
            wbar0 = 0.66567488;  %0.49949686;  %theta3: 0.66567488; 
            chi0_min = 0.8;
            chi0_max = 1.1;
            wbar0_min = 0.4;
            wbar0_max = 0.6;
    case 4  %Woodford 
            switch mymoment
                        case 'hazard'
                            menucost0 = 0.096444681;  %0.23903585; %0.04;  %kurt/Theta2: 0.041741066; %hazard: 1.1723458;  %kurt: 1.2040397; 
                            std_eA0 = 0.03898251;  %0.055176423; %0.058854404; %kurt/Theta2: 0.0500062;  %hazard: 0.056698077;  %kurt: 0.056729418;
                            alpha0  = 0.28239298;  %0.69354316;  %%5; %kurt/Theta2: 0.028145528;%hazard: 6.6177836;  %kurt: 0.20718998;
                            chi0 = 0.67800877;  %0.90292731;  %%0.68601033; %kurt/Theta2: 0.68858532; %hazard: 0.99722296;  %kurt: 0.99339566;
                            wbar0 = 0.6621076;  %0.41774662;  %0.31574312; %kurt/Theta2: 0.32700956;  %%hazard: 0.43809254;  %kurt: 0.43568317;
                            std_m0 = 0.0013480617; %0.0033991633; %monthly inflation matched: 0.0027272733;
                            std_A0 = 0.0094349537;
                        case 'kurtosis'
                            menucost0 = 0.031807171; %%hazard/Theta2: 0.23903585; %0.04;  %kurt/Theta2: 0.041741066; %hazard: 1.1723458;  %kurt: 1.2040397; 
                            std_eA0 = 0.036835504;  %hazard/Theta2: 0.055176423; %0.058854404; %kurt/Theta2: 0.0500062;  %hazard: 0.056698077;  %kurt: 0.056729418;
                            alpha0  = 0.023255054; %hazard/Theta2: 0.69354316;  %%5; %kurt/Theta2: 0.028145528;%hazard: 6.6177836;  %kurt: 0.20718998;
                            chi0 = 0.43121035; %hazard/Theta2: 0.90292731;  %%0.68601033; %kurt/Theta2: 0.68858532; %hazard: 0.99722296;  %kurt: 0.99339566;
                            wbar0 = 0.42738617;  %hazard/Theta2: 0.41774662;  %0.31574312; %kurt/Theta2: 0.32700956;  %%hazard: 0.43809254;  %kurt: 0.43568317;
                        case 'stdpcdist'
                            menucost0 = 0.098887803;  %0.10233266; %%hazard/Theta2: 0.23903585; %0.04;  %kurt/Theta2: 0.041741066; %hazard: 1.1723458;  %kurt: 1.2040397; 
                            std_eA0 = 0.050141763;  %0.050371213;  %hazard/Theta2: 0.055176423; %0.058854404; %kurt/Theta2: 0.0500062;  %hazard: 0.056698077;  %kurt: 0.056729418;
                            alpha0  = 0.073649332;  %0.080425336; %hazard/Theta2: 0.69354316;  %%5; %kurt/Theta2: 0.028145528;%hazard: 6.6177836;  %kurt: 0.20718998;
                            chi0 = 0.9074341;  %0.90736477; %hazard/Theta2: 0.90292731;  %%0.68601033; %kurt/Theta2: 0.68858532; %hazard: 0.99722296;  %kurt: 0.99339566;
                            wbar0 = 0.42853758;  %0.4282449;  %hazard/Theta2: 0.41774662;  %0.31574312; %kurt/Theta2: 0.32700956;  %%hazard: 0.43809254;  %kurt: 0.43568317;
            end
            chi0_min = 0.85;
            chi0_max = 0.95;
            wbar0_min = 0.40;
            wbar0_max = 0.45;
    case 5  %linear hazard
            switch mymoment
                        case 'hazard'
                            std_eA0 = 0.037482317;  
                            alpha0  = 0.71650986; 
                            chi0 = 0.67646765; 
                            wbar0 = 0.6631284;
                            omega0 = 0.06113182; 
                            chi0_min = 0.85;
                            chi0_max = 0.95;
                            wbar0_min = 0.40;
                            wbar0_max = 0.45;
                        case 'kurtosis'
                            std_eA0 = 0.035485222;  
                            alpha0  = 1.3725519; 
                            chi0 = 0.43113057; 
                            wbar0 = 0.42763255;
                            omega0 = 0.02859706; 
                            chi0_min = 0.85;
                            chi0_max = 0.95;
                            wbar0_min = 0.40;
                            wbar0_max = 0.45;
            end
    case 6  %asymmetric linear hazard
                    switch switch_distribution
                        case 'gaussian'
                            std_eA0 = 0.037802453;  
                            alpha0  = 0.8183596; 
                            alphapos0  = 0.48743093; 
                            chi0 = 0.67673962; 
                            wbar0 = 0.66300079;
                            omega0 = 0.065437287;
                        case 'laplace' %Not done
                            std_eA0 = 0.053879927 ; %0.056373756;  
                            alpha0  = 0.8183585; %0.59430218; 
                            alphapos0  = 0.48743507; %0.59430218; 
                            chi0 = 0.89899678; %0.90451547; 
                            wbar0 = 0.42340555; %0.41965358;
                            omega0 = 0.066476407; %0.086127045;
                        case 'mixed' %Not done
                            std_eA0 = 0.081656119; %0.056373756;  
                            alpha0  = 0.75642074; %0.59430218; 
                            alphapos0  = 0.38390946; %0.59430218; 
                            chi0 = 0.89628242; %0.90451547; 
                            wbar0 = 0.41588859; %0.41965358;
                            omega0 = 0.090739527; %0.086127045;
                            lam0 = 1;
                    end
            std_m0 = 0.001333697; %monthly inflation matched: 0.0027272733;  
            std_A0  = 0.0084540718;
            chi0_min = 0.85;
            chi0_max = 0.95;
            wbar0_min = 0.40;
            wbar0_max = 0.45;
    case 7  %uniform random menu cost
            std_eA0 = 0.038525554;  
            alpha0  = 0.59391691; 
            chi0 = 0.67794289; 
            wbar0 = 0.66309512;
            omega0 = 0.0818877; 
            chi0_min = 0.85;
            chi0_max = 0.95;
            wbar0_min = 0.40;
            wbar0_max = 0.45;
end

mu_ = 0; % annual net inflation rate

ksi=NaN; %#ok<*NASGU> % initialize legacy code params

mfreq = 12; % model frequency: 4 quarterly; 12 monthly; 365/7 weekly
   
% Macro parameters
beta  = (1.03)^(-1/mfreq); % time discount factor
gam  = 2;  %2;                 % CRRA coefficient
theta = 3;  %3;                 % Dixit-Stiglitz elast. of subst.
nu    = 1;                 % log money demand parameter
eta   = 0/3; %1/3;               % degree of real rigidity (decreasing returns)

% Autocorrelations 
rho  = 1;   % Idiosyncratic productivity process  

mu = (1 + mu_/100).^(1/mfreq); % gross inflation at model frequency

jacstep = sqrt(eps); 
par_.jacstep = jacstep;
par_.switch_shock = switch_shock;
par_.switch_distribution = switch_distribution;